Modèle linéaire multiple

2eme-FA-EMS - BUT SD - E. Anakok

0.1 Ce qu’il faut retenir de ce cours

Modélisation probabiliste du modèle linéaire

\(y_i\) est la réalisation d’une variable aléatoire \(Y_i\) telle que pour \(1 \leq i \leq n\):

\[Y_i = \beta + \sum_{j=1}^p\alpha_j x_{i,j} + E_i, \quad E_i\overset{i.i.d.}\sim{\cal N}(0,\sigma^2) \]

Écriture matricielle du modèle

Le modèle linéaire multiple peut s’écrire matriciellement

\[ Y = X \theta + E \]

Estimateur des moindres carrés de \(\theta\)

\[\widehat{\theta} =\underset{\theta\in\mathbb{R}^{p+1}}{argmin} ||Y-X\theta||^2\]

Si \((X^t\,X)\) est inversible, alors

\[\textcolor{red}{\hat{\theta}=(X^t\,X)^{-1}\,X^t\,Y}\]

Vérification des hypothèses : Graphe de diagnostic, corrélation et vif

Warning: `fortify(<lm>)` was deprecated in ggplot2 4.0.0.
ℹ Please use `broom::augment(<lm>)` instead.
ℹ The deprecated feature was likely used in the ggfortify package.
  Please report the issue at <https://github.com/sinhrks/ggfortify/issues>.

Warning: `aes_string()` was deprecated in ggplot2 3.0.0.
ℹ Please use tidy evaluation idioms with `aes()`.
ℹ See also `vignette("ggplot2-in-packages")` for more information.
ℹ The deprecated feature was likely used in the ggfortify package.
  Please report the issue at <https://github.com/sinhrks/ggfortify/issues>.

Warning: Using `size` aesthetic for lines was deprecated in ggplot2 3.4.0.
ℹ Please use `linewidth` instead.
ℹ The deprecated feature was likely used in the ggfortify package.
  Please report the issue at <https://github.com/sinhrks/ggfortify/issues>.

\[ VIF(x_i) = \frac{1}{1-R^2(x_i)} \]

0.2 Problématique biologique

Données On a pour 20 brèmes péchées dans le lac Laengelmavesi en Finland leurs poids (en gramme) et leurs longeurs (en cm), leurs largeurs (en cm) leurs épaisseur (en cm).

Objectif

Peut-on prédire le poids des poissons par leurs largeurs, leurs longueurs et leurs épaisseurs ?

Tester le caractère significatif de la liaison : Y’a t’il un lien significatif entre le poids des poissons et leurs largeurs, leurs longueurs et leurs épaisseurs ?

Y’a t’il un lien entre le poids du poisson et sa largeur après avoir pris en compte sa longueur et son épaisseur ?

1 Écriture et remarque sur le modèle

1.1 Écriture du modèle dans cet exemple

1.1.0.1 Notations

On a \(n = 20\) observations. On note, pour \(1 \leq i \leq 20\)

\(x_{1,i}\) la mesure de la longueur du poisson \(i\).
\(x_{2,i}\) la mesure de la largeur du poisson \(i\).
\(x_{3,i}\) la mesure de l’épaisseur du poisson \(i\).
\(y_i\) la mesure du poids du poisson \(i\).

Définition : Modèle de régression linéaire multiple

On suppose que \(y_i\) est la réalisation d’une variable aléatoire \(Y_i\) telle que pour \(1\leq i \leq n\)

\[Y_i = \alpha_1 x_{1,i} + \alpha_2 x_{2,i} + \alpha_3 x_{3,i} + \beta + E_i\]

\(\alpha_1\) (resp \(\alpha_2\) , \(\alpha_3\) ) est un paramètre inconnu, l’effet de la longueur (resp largeur, épaisseur) sur le poids;
\(\beta\) est un paramètre inconnu;
\(E_i\) une variable aléatoire, \(E_i \overset{i.i.d.}\sim \mathcal{N}(0, \sigma^2)\)
Les \(x_i\) sont linéairement indépendants

1.2 Écriture plus générale du modèle

1.2.0.1 Notations

On a \(n = 20\) observations. On note, pour \(1 \leq i \leq 20\)

\(x_{1}, \dots, x_{p}\), \(p\) variables explicatives et \(y\) une variable réponse
\(y_i\) la mesure de la variable réponse de l’individu \(i\).
\(\forall j \in \{ 1,\dots, p \}\), \(x_{j,i}\) la valeur de la variable \(x_j\) de l’individu \(i\)

Définition : Modèle de régression linéaire multiple

On suppose que \(y_i\) est la réalisation d’une variable aléatoire \(Y_i\) telle que pour \(1 \leq i \leq n\) : \[Y_i = \sum_{j=1}^p \alpha_j x_{j,i} + \beta + E_i \]

\(\alpha_j\) est un paramètre inconnu, l’effet de la variable \(x_j\) sur le y.
\(\beta\) est un paramètre inconnu.
\(E_i\) une variable aléatoire, \(E_i \overset{i.i.d.}\sim \mathcal{N}(0, \sigma^2)\).
Les \(x_i\) sont linéairement indépendants.

1.3 Deux écritures équivalentes

Remarques

Les deux modèles suivant sont équivalent :

\(\forall i \in \{1, \dots, n\} \;\;\)

\[ Y_i = \sum_{j=1}^p \alpha_j x_{j,i} + \beta + E_i,\quad \ E_i \overset{i.i.d.}\sim \mathcal{N}(0, \sigma^2)\]

\(\forall i \in \{1, \dots, n\} \;\;\)

\[ Y_i \overset{i.i.d.}\sim \mathcal{N}\left(\sum_{j=1}^p \alpha_j x_{j,i} + \beta, \sigma^2\right)\]

1.4 Dans l’exemple des poissons

Remarques

Les deux modèles suivant sont équivalent :

\(\forall i \in \{1, \dots, n\} \;\;\) \[ Y_i = \alpha_1 x_{1,i} + \alpha_2 x_{2,i} + \alpha_3 x_{3,i} + \beta + E_i\]

\[E_i \overset{i.i.d.}\sim \mathcal{N}(0, \sigma^2)\]

\(\forall i \in \{1, \dots, n\} \;\;\) \[Y_i \overset{i.i.d.}\sim \mathcal{N}\left(\alpha_1 x_{1,i} + \alpha_2 x_{2,i} + \alpha_3 x_{3,i} + \beta, \sigma^2\right)\]

1.5 Remarques sur les paramètres

Remarques

Dans le modèle \(\forall i \in \{1, \dots, n\} \;\;\)

\[ Y_i = \sum_{j=1}^p \alpha_j x_{j,i} + \beta + E_i,\quad\ E_i \overset{i.i.d.}\sim \mathcal{N}(0, \sigma^2)\]

\(\alpha_j\) représente l’accroissement de \(Y_i\) correspondant à l’accroissement d’une unité sur \(x_j\) si les autres variables explicatives sont fixées. Dans l’exemple des poissons : \(\alpha_1\) represente l’accroissement du poids correspondant à l’accroissement d’1 cm sur la longeurs du poisson si la largeur et l’épaisseur son fixé.

Ce modèle est de dimension \(p + 1\). (Avec \(p+1\) paramètres d’espérance à estimer \(p\) pour \(\alpha\) et \(1\) pour\(\beta\))

\(\forall i \in \{1, \dots, n\}\) \(Y_i\) se décompose en \(\color{red}{Y_i} = {\color{blue}{\underbrace{\displaystyle\sum_{j=1}^p \alpha_j x_{j,i} + \beta}_{{déterministe}{}}}} + \color{red}{\overbrace{E_i}^{aléatoire}}\)

2 Écriture matricielle du modèle

2.1 Rappel sur les matrices

Qu’est ce qu’une matrice ?
La multiplication matricielle
La transposée d’une matrice
L’inverse d’une matrice
Petits topo sur les vecteurs

2.2 Rappels d’algèbre

Rappel : Matrice

\(A = (a_{i,j})_{1\leq i \leq p, 1\leq j \leq q}\) est une matrice à coefficients réels de format \((p \times q)\)

\[A = \begin{bmatrix} a_{11} & a_{12} & \ldots & a_{1q}\\ a_{21} & a_{22} & \ldots & a_{2q}\\ \vdots & \vdots & & \vdots \\ a_{p1} & a_{p2} & \ldots & a_{pq} \end{bmatrix}\]

Rappel : Transposée d’une matrice

Si \(A = (a_{i,j})_{1\leq i \leq p, 1\leq j \leq q}\) est une matrice \((p \times q)\), on note alors \(A^t\) sa transposée, la matrice \((q \times p)\) définie par

\[A^t = \begin{bmatrix} a_{11} & a_{21} & \ldots & a_{p1}\\ a_{12} & a_{22} & \ldots & a_{p2}\\ \vdots & \vdots & & \vdots \\ a_{1q} & a_{2q} & \ldots & a_{pq} \end{bmatrix}\]

Rappel : Matrice symétrique

\(A\) est une matrice symétrique si :

\(A\) est une matrice carré (\(p=q\))
\(A^t = A\)

Rappel : Transposée du produit

Soit \(A\) et \(B\) deux matrices. On a

\[ (A B)^t = B^t A^t\]

2.3 Rappels d’algèbre : Vecteur

Rappel : Vecteurs

On note V un vecteur à \(p\) éléments par

\[V = \begin{bmatrix} v_{1} \\ v_{2}\\ \vdots \\ v_{p} \end{bmatrix}= \left[{v_1, v_2, \ldots, v_p}\right]^t \]

2.4 Rappels d’algèbre : Vecteur

Rappel : Produit matrice-vecteur

Si \(V\) est un vecteur à \(q\) éléments et \(A\) une matrice \((p \times q)\) alors, \(U=AV\) est un vecteur à \(p\) éléments

\[\begin{align} U &= \begin{bmatrix} u_{1} \\ u_{2}\\ \vdots \\ u_{p} \end{bmatrix} = AV \\ \\ AV&=\begin{bmatrix} a_{11} & a_{12} & \ldots & a_{1q}\\ a_{21} & a_{22} & \ldots & a_{2q}\\ \vdots & \vdots & & \vdots \\ a_{p1} & a_{p2} & \ldots & a_{pq} \end{bmatrix}\begin{bmatrix} v_{1} \\ v_{2}\\ \vdots \\ v_{p} \end{bmatrix}\\ \\AV&= \begin{bmatrix} \sum_{j=1}^{q} a_{1j} v_j \\ \sum_{j=1}^{q} a_{2j} v_j\\ \vdots \\ \sum_{j=1}^{q} a_{pj} \end{bmatrix} \end{align}\]

2.5 Rappels d’algèbre : Identité et inverse

Rappel : matrice identité

La matrice identité \(I_n\) est la matrice carré à \(n\) lignes et \(n\) colonnes, dont les coefficients diagonaux valent \(1\) et les autres valent \(0\).

\[I_n = \begin{bmatrix} 1 & 0 & \ldots & 0\\ 0 & 1 & \ldots & 0\\ \vdots & \vdots & & \vdots \\ 0 & 0 & \ldots & 1 \end{bmatrix}\]

Rappel : Inverse d’une matrice

Si \(A = (a_{i,j})_{1\leq i \leq p, 1\leq j \leq p}\) est une matrice carrée \((p \times p)\) inversible, alors on note son inverse \(A^{-1}\) la matrice \((p \times p)\) telle que

\[A\,A^{-1} = A^{-1}\,A = I_p \]

Rappel : Transposée de l’inverse

Si \(A\) est une matrice inversible, alors

\[(A^{-1})^t=(A^t)^{-1}\]

2.6 Vecteurs aléatoires

Rappel : Vecteurs aléatoires

Soient \(V_1,V_2,\dots,V_n\), \(n\) variables aléatoires réelles, alors

\[V = \begin{bmatrix} V_1 \\ V_2 \\ \vdots \\ V_n \end{bmatrix} \]

Rappel : Espérance d’un vecteur aléatoire

\(V\) a pour espérance

\[\mathbb{E}[V] = \begin{bmatrix} \mathbb{E}[V_1] \\ \mathbb{E}[V_2] \\ \vdots \\ \mathbb{E}[V_n] \end{bmatrix}\]

Rappel : Matrice de variance-covariance

\(V\) a pour matrice de variance-covariance (symétrique)

\[\mathbb{V}[V] = \mathbb{E}\left[\left(V- \mathbb{E}[V] \right)\left(V- \mathbb{E}[V] \right)^t \right] =\]

\[\begin{bmatrix} \mathbb{V}[V_1] & cov(V_1,V_2) & \dots& Cov(V_1,V_n) \\ cov(V_2,V_1) & \mathbb{V}[V_2] & \dots& Cov(V_2,V_n) \\ \vdots & \vdots & \ddots & \vdots \\ cov(V_n,V_1) & cov(V_n,V_2) & \dots& \mathbb{V}[V_n]\end{bmatrix}\]

2.7 Propriétés

Théorème

Soient \(\textcolor{blue}{A}\) et \(\textcolor{blue}{B}\) des matrices à coefficients non aléatoires et soit \(\textcolor{red}{V}\) un vecteur aléatoire, alors

\[\begin{align} \mathbb{E}[\textcolor{blue}{A}\textcolor{red}{V}\textcolor{blue}{B}] &= \textcolor{blue}{A} \mathbb{E}[\textcolor{red}{V}]\textcolor{blue}{B}\\ \mathbb{V}[\textcolor{blue}{A}\textcolor{red}{V}] & = \textcolor{blue}{A} \mathbb{V}[\textcolor{red}{V}] \textcolor{blue}{A^t} \end{align}\]

Rappel : Vecteurs gaussiens

Si \(V_1, V_2, \ldots, V_n\) sont gaussiennes et indépendantes, alors \(V = [V_1, \ldots, V_n]^t\) est un vecteur gaussien.

Si \(V_1, V_2, \ldots, V_n \overset{i.i.d}\sim \mathcal{N}(0,\sigma ^2)\), alors

\[V \sim \mathcal{N}\left(\vec{0}_n,\,\sigma^2 I_n\right)\] où \(I_n\) est la matrice identité d’ordre \(n\).

2.8 Écriture matricielle du modèle

Théorème

\[\left\{\begin{matrix} Y_1 = \beta\ +\ \alpha_1\,x_{1,1}\ +\ \alpha_2\,x_{2,1}\ +\ \ldots\ +\ \alpha_{p}\,x_{p, 1}\ +\ E_1\\ Y_2 = \beta\ +\ \alpha_1\,x_{1, 2}\ +\ \alpha_2\,x_{2,2}\ +\ \ldots\ +\ \alpha_{p}\,x_{p, 2}\ +\ E_2\\ \vdots \\ Y_n = \beta\ +\ \alpha_1\,x_{1,n}\ +\ \alpha_2\,x_{2, n}\ +\ \ldots\ +\ \alpha_{p}\,x_{p,n}\ +\ E_n \end{matrix} \right.\]

peut se réécrire sous la forme

\[ Y = X \theta + E \]

avec

\[Y =\begin{bmatrix} Y_1\\ Y_2\\ \vdots \\ Y_n\end{bmatrix},\quad X = \begin{bmatrix} 1 & x_{1,1} & \ldots & x_{p, 1}\\ 1 & x_{1, 2} & \ldots & x_{p, 2}\\ \vdots & \vdots & & \vdots \\ 1 & x_{1, n} & \ldots & x_{p, n} \end{bmatrix}\]

\[\theta = \begin{bmatrix} \beta\\ \alpha_1\\ \alpha_2\\ \vdots \\ \alpha_{p}\end{bmatrix},\quad E = \begin{bmatrix}E_1\\ E_2\\ \vdots \\ E_n\end{bmatrix}\]

2.9 Écriture matricielle du modèle

Propriété

\(X\) est de taille \(n\times (p+1)\)

\(\theta\) est de dimension \((p+1)\times 1\) : vecteur des paramètres d’espérance

\(E\) est de dimension \((n\times 1)\) : est un vecteur Gaussien tel que

\[E \sim \mathcal{N}\left(\vec{0}_n,\,\sigma^2 I_n\right) \]

\(Y\) est de dimension \((n\times 1)\) : est un vecteur Gaussien tel que

\[ Y \sim \mathcal{N}\left(X\theta,\,\sigma^2 I_n\right)\]

Remarques

Si les vecteurs colonnes de \(X\) sont linéairement indépendants, \(X\) est de rang \(p+1\) et \(X^t X\) est inversible.

2.10 Exercice : écrire le modèle linéaire univarié sous forme matricielle

Le modèle de régression simple \[Y_i=\alpha x_i+\beta+E_i\] s’écrit matriciellement \(Y=X\theta+E\) avec \[Y = \begin{bmatrix} Y_1\\ Y_2\\ \vdots \\ Y_n\end{bmatrix},\quad X = \begin{bmatrix} 1 & x_{1} \\ 1 & x_{2} \\ \vdots & \\ 1 & x_{n} \\ \end{bmatrix}\]

\[\theta\ =\ \begin{bmatrix}\beta\\ \alpha\end{bmatrix},\quad E\ =\ \begin{bmatrix} E_1\\ E_2\\ \vdots \\ E_n\end{bmatrix}\]

2.11 Écriture du modèle linéaire simple

Propriété

\(X\) est de taille \((n\times 2)\)
\(\theta\) est de dimension \((2 \times 1)\) : vecteur des paramètres d’espérance
\(E\) est de dimension \((n\times 1)\) : est un vecteur Gaussien tel que

\[E \sim \mathcal{N}\left(\vec{0}_n,\,\sigma^2 I_n\right) \]

\(Y\) est de dimension \((n\times 1)\) : est un vecteur Gaussien tel que

\[ Y \sim \mathcal{N}\left(X\theta,\,\sigma^2 I_n\right)\]

Remarques

Si les \(x_i\) ne sont pas tous égaux, les deux vecteurs colonnes de \(X\) sont linéairement indépendants et \(X\) est de rang \(p = 2\).

2.12 Exemple sur les poissons

2.12.0.1 Notations

On a \(n = 20\) observations. On note, pour \(1 \leq i \leq 20\)

\(x_{1,i}\) la mesure de la longueur du poisson \(i\).
\(x_{2,i}\) la mesure de la largeur du poisson \(i\).
\(x_{3,i}\) la mesure de l’épaisseur du poisson \(i\).
\(y_i\) la mesure du poids du poisson \(i\).

Définition : Modèle de régression linéaire multiple

On suppose que \(y_i\) est la réalisation d’une variable aléatoire \(Y_i\) telle que pour \(1 \leq i \leq n\) :

\[Y_i = \alpha_1 x_{1,i} + \alpha_2 x_{2,i} + \alpha_3 x_{3,i} + \beta + E_i,\quad \]

\(\alpha_1\) (resp \(\alpha_2\) , \(\alpha_3\) ) est un paramètre inconnu, l’effet de la longueur (resp largeur, épaisseur) sur le poids;
\(\beta\) est un paramètre inconnu;
\(E_i\) une variable aléatoire, \(E_i \overset{i.i.d.}\sim \mathcal{N}(0, \sigma^2)\)
Les \(x_i\) sont linéairement indépendants

2.13 Comment écrire le modèle :

\[\left\{\begin{matrix} Y_1 = \beta\ +\ \alpha_1\,x_{1,1}\ +\ \alpha_2\,x_{2,1}\ +\ \alpha_{3}\,x_{3, 1}\ +\ E_1\\ Y_2 = \beta\ +\ \alpha_1\,x_{1, 2}\ +\ \alpha_2\,x_{2,2}\ +\ \alpha_{3}\,x_{3, 2}\ +\ E_2\\ \vdots \\ Y_n = \beta\ +\ \alpha_1\,x_{1,n}\ +\ \alpha_2\,x_{2, n}\ +\ \alpha_{3}\,x_{3,n}\ +\ E_n \end{matrix} \right.\]

peut se réécrire sous la forme

\[ Y = X \theta + E \]

avec

\[Y =\begin{bmatrix} Y_1\\ Y_2\\ \vdots \\ Y_n\end{bmatrix},\quad X = \begin{bmatrix} 1 & x_{1,1} & x_{2,1} & x_{3, 1}\\ 1 & x_{1, 2} & x_{2,1} & x_{3, 2}\\ \vdots & \vdots & & \vdots \\ 1 & x_{1, n} & x_{2,1} & x_{3, n} \end{bmatrix}\]

\[\theta = \begin{bmatrix} \beta\\ \alpha_1\\ \alpha_2\\\vdots\\ \alpha_{3}\end{bmatrix},\quad E = \begin{bmatrix}E_1\\ E_2\\ \vdots \\ E_n\end{bmatrix}\]

2.14 Avec les données :

\[X = \]

Intercept	Length1	Height	Width
1	23.2	11.5200	4.0200
1	24.0	12.4800	4.3056
1	23.9	12.3778	4.6961
1	26.3	12.7300	4.4555
1	26.5	12.4440	5.1340
1	26.8	13.6024	4.9274
1	26.8	14.1795	5.2785
1	27.6	12.6700	4.6900
1	27.6	14.0049	4.8438
1	28.5	14.2266	4.9594

3 Estimation des paramètres d’espérance du modèle

3.1 Estimateurs des moindres carrés

Objectif

On note \(\theta= [\beta,\alpha_1,...,\alpha_p]^t\).

Comme dans la regression linéaire simple on cherche à minimiser :

\[J(\theta) = \sum_{i = 1}^n \left( y_i -(\beta + \sum_{j=1}^p \alpha _j x_{j,i}) \right)^2\] Écriture matricielle ?

3.2 Rappels d’algèbre : norme

Rappel : norme

Soit \[V\ =\ \begin{bmatrix} v_{1} \\ v_{2}\\ \vdots \\ v_{n}\end{bmatrix} \ =\ [v_1, \ldots, v_n]^t\] un vecteur de \(\mathbb{R}^n\). On appelle norme de \(V\) le réel

\[\displaystyle V^t V \ =\ \sum_{i=1}^{n} v_i^2\ =\ ||V||^2 \]

Remarque

\[V^t V \neq V V^t\]

\[V V^t = \begin{bmatrix} v_1^2 & v_1\,v_2 & \ldots & v_1\,v_n\\ v_2\,v_1 & v_2^2 & \ldots & v_2\,v_n\\ \vdots & \vdots & & \vdots \\ v_n\,v_1 & v_n\,v_2 & \ldots & v_n^2 \end{bmatrix}\]

3.3 Estimateurs des moindres carrés

Objectif

On souhaite trouver \(\theta= (\beta,\alpha_1,...\alpha_p)^t\) qui minimise

\[\begin{align}J(\theta)&=\sum_{i = 1}^n \left( y_i -(\beta + \sum_{j=1}^p \alpha _j x_{j,i}) \right)^2 \\&= (Y-X\theta)^t(Y-X\theta) \\&= ||Y-X\theta||^2\end{align}\]

Théorème

Le \(\theta\) optimal est la solution du système de \(p+1\) équations à \(p+1\) inconnues

\[ \left\{\begin{array}{ccc} \frac{\partial \sum_{i = 1}^n \left( y_i -(\beta + \sum_{j=1}^p \alpha _j x_{j,i}) \right)^2 }{\partial \beta} & = & 0 \\ % \frac{\partial \sum_{i = 1}^n \left( y_i -(\beta + \sum_{j=1}^p \alpha _j x_{j,i}) \right)^2 }{\partial \alpha_1} & = & 0\\ \vdots & & \\ \frac{\partial \sum_{i = 1}^n \left( y_i -(\beta + \sum_{j=1}^p \alpha _j x_{j,i}) \right)^2 }{\partial \alpha_p} & = & 0 \end{array}\right. \]

3.4 À la recherche de \(\widehat{\theta}\)

Estimateur des moindres carrés

\[\begin{align}\widehat{\theta} &= \underset{\theta\in\mathbb{R}^{p+1}}{argmin} (Y-X\theta)^t(Y-X\theta) \\ \widehat{\theta}&=\underset{\theta\in\mathbb{R}^{p+1}}{argmin} ||Y-X\theta||^2\end{align}\]

\(\widehat{\theta} = [B,A_1, \dots, A_p]^t\) où \(B\) est l’estimateur de \(\beta\) et \(A_j\) l’estimateur de \(\alpha_j\)

Écriture matricielle des équations

\[\begin{align} J(\theta) &= (Y-X\theta)^t(Y-X\theta)\\ &=(Y^t-\theta^tX^t)(Y-X\theta)\\ &= Y^tY - Y^tX\theta -\theta^tX^t Y +\theta^tX^tX\theta \\ &= Y^tY - 2 Y^tX\theta +\theta^tX^tX\theta \end{align}\]

En dérivant par rapport à \(\theta\) on obtient : \[-2Y^tX + 2\theta^tX^tX \]

Le système \(-2Y^tX + 2\theta^tX^tX=0\) est bien le même que les \(p+1\) équations de la slide précédente

3.5 À la recherche de \(\widehat{\theta}\)

Remarques

\(X\) est de rang \(p+1\) donc \(X^tX\) est inversible.

Fin du calcul

\[\begin{align} -2Y^tX + 2\widehat{\theta}^tX^tX = 0 &\iff \widehat{\theta}^tX^tX = Y^tX \\ \iff X^tX \widehat{\theta} = X^tY & \iff \widehat{\theta} = (X^tX)^{-1} X^tY \end{align}\]

Théorème

\[\widehat{\theta} =\underset{\theta\in\mathbb{R}^{p+1}}{argmin} ||Y-X\theta||^2\] a pour solution

\[\widehat{\theta} = (X^tX)^{-1} X^tY\]

3.6 Estimateurs des paramètres d’espérance

Objectif : Estimateur des moindres carrés de \(\theta\)

\[\widehat{\theta} =\underset{\theta\in\mathbb{R}^{p+1}}{argmin} ||Y-X\theta||^2\]

Théorème

Si \(\hat{\theta}= [B, A_1, \ldots, A_{p}]^t\) est solution du système et que \((X^t\,X)\) est inversible, alors

\[\hat{\theta}=(X^t\,X)^{-1}\,X^t\,Y\]

Paramètres	Estimateurs	Estimations
\([\beta,\alpha_1,...\alpha_p]\)	\([B, A_1, \ldots, A_{p}]\)	\([b,a_1,\cdots,a_p]\)

4 Validation des hypothèses et estimation avec R

4.1 Les 4 graphiques

4.2 Multicolinéarité des \(x_i\)

On vérifie la corrélation entre les \(x_i\)

library(corrplot)
M = cor(fish[c("Length1","Height","Width","Weight")])
corrplot.mixed(M)

4.3 Vérification avec le VIF

Définition : VIF

\[ VIF(x_i) = \frac{1}{1-R^2(x_i)} \] où \(R^2(x_i)\) est le \(R^2\) du modèle linéaire dans lequel on explique \(x_i\) par toutes les autres covariables (\(x_j\))

4.4 Exemple du VIF sur les poissons :

Calcul ‘à la main’

mod_length<- lm(Length1 ~ Height + Width , data =  fish)
vif_length <- 1/ ( 1-summary(mod_length)$r.squared)
print(vif_length)

[1] 8.952978

Avec la fonction VIF :

library(car)
vif(mod)

  Length1    Height     Width 
 8.952978 12.123683  7.451746

En pratique : on retire les covariables qui ont un VIF supérieur à 10 (ou 5 si on veut être plus strict)

4.5 Estimations des paramètres sur notre exemple :

mod <- lm(Weight ~Length1 + Height + Width , data =  fish)
summary(mod)


Call:
lm(formula = Weight ~ Length1 + Height + Width, data = fish)

Residuals:
     Min       1Q   Median       3Q      Max 
-183.872  -17.044   -0.495   24.591  101.032 

Coefficients:
            Estimate Std. Error t value Pr(>|t|)    
(Intercept) -1048.57     187.86  -5.582 4.13e-05 ***
Length1        18.96      13.26   1.430    0.172    
Height         48.50      28.31   1.713    0.106    
Width          66.72      51.73   1.290    0.215    
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Residual standard error: 63.03 on 16 degrees of freedom
Multiple R-squared:  0.808, Adjusted R-squared:  0.772 
F-statistic: 22.44 on 3 and 16 DF,  p-value: 5.624e-06

4.6 Sur notre exemple (sans la fonction lm)

On récupère les matrices \(X\) et \(Y\) :

X1 <- fish[, c("Length1",  "Height",  "Width")] %>% 
  as.matrix() 
X <- cbind(Intercept = 1,X1)
Y <- fish[,"Weight"] %>% 
  as.matrix()

4.7 Sur notre exemple (sans la fonction lm)

On cherche maintenant \(\widehat{\theta}=(X^t\,X)^{-1}\,X^t\,Y\).

Le produit matriciel \(AB\) s’obtient en faisant A %*% B, La transposée d’une matrice A s’obtient avec la fonction t(A) et son inverse avec la fonction solve(A).

\(X^t\,X\) s’obtient donc :

t(X) %*% X

\((X^t\,X)^{-1}\) s’obtient donc :

solve(t(X) %*% X)

On a finalement \(\widehat{\theta}\) :

solve(t(X) %*% X)%*% t(X)%*% Y

4.8 On retrouve bien les mêmes valeurs

coefficients(mod)

(Intercept)     Length1      Height       Width 
-1048.57013    18.95867    48.49646    66.71559

solve(t(X) %*% X)%*% t(X)%*% Y

                 [,1]
Intercept -1048.57013
Length1      18.95867
Height       48.49646
Width        66.71559

4.9 Ce qu’il faut retenir de ce cours

Modélisation probabiliste du modèle linéaire

\(y_i\) est la réalisation d’une variable aléatoire \(Y_i\) telle que pour \(1 \leq i \leq n\):

\[Y_i = \beta + \sum_{j=1}^p\alpha_j x_{i,j} + E_i, \quad E_i\overset{i.i.d.}\sim{\cal N}(0,\sigma^2) \]

Écriture matricielle du modèle

Le modèle linéaire multiple peut s’écrire matriciellement

\[\textcolor{red}{Y = X \theta + E}\]

Estimateur des moindres carrés de \(\theta\)

\[\widehat{\theta} =\underset{\theta\in\mathbb{R}^{p+1}}{argmin} ||Y-X\theta||^2\]

Si \((X^t\,X)\) est inversible, alors \(\hat{\theta}=(X^t\,X)^{-1}\,X^t\,Y\)

Vérification des hypothèses : Graphe de diagnostic, corrélation et vif

\[ VIF(x_i) = \frac{1}{1-R^2(x_i)} \]

5 Estimateur de la variance

5.1 Ce qu’il faut retenir de ce cours 1/2

Estimateur de \(\sigma^2\) et loi

\[{S^2}_{(M_{p+1})}=\frac{SCR(M_{p+1})}{n-(p+1)}=\frac{||Y - X\,\widehat{\theta}||^2}{n-p-1}\]

\[\frac{(n-p-1){S^2}_{(M_{p+1})}}{\sigma^2}\sim\chi^2(n-p-1)\]

Loi de \(\widehat{\theta}\) et \(\widehat{\theta}_k\)

\[\widehat{\theta} \sim \mathcal{N}\left(\theta,\sigma^2(X^tX)^{-1}\right)\] \[\widehat{\theta}_k \sim \mathcal{N}\left(\theta_k,\sigma^2c_{kk}\right)\text{ avec } c_{kk} = (X^tX)^{-1}_{k,k}\]

\[\frac{\widehat{\theta}_k- \theta_k}{\sqrt{S^2_{(M_{p+1})}c_{kk}}}\sim\mathcal{T}(n-p-1) \]

Intervalle de confiance de \(\theta_k\)

\[\begin{align} IC_{1-\delta}(\theta_k) = &\left[\hat{\theta}_k-t_{1-\frac{\delta}{2}} \sqrt{{S^2}_{(M_{p+1})}c_{kk}}\right.;\\ &\;\left.\hat{\theta}_k+t_{1-\frac{\delta}{2}} \sqrt{{S^2}_{(M_{p+1})}c_{kk}}\;\right]\end{align}\]

5.2 Ce qu’il faut retenir de ce cours 2/2

Test de Student de la nullité de la pente

\[H_0:\alpha_k=0\]

\[H_1:\alpha_k\neq 0\]

Test de Fisher de comparaison de \(M_1\) contre \(M_{p+1}\)

\[H_0\; :\; \text{modèle}\; M_1:\] \[Y_i=\beta+E_i,\quad E_i\; \overset{i.i.d.}{\sim} {\cal N}(0,\sigma^2)\]

\[H_1\; :\; \text{modèle}\; M_{p+1}:\] \[Y_i= \beta + \sum_{j=1}^p\alpha_j x_{j,i}+E_i,\; E_i\; \overset{i.i.d.}{\sim} {\cal N}(0,\sigma^2)\]

Test de Fisher de comparaison de \(M_{q+1}\) contre \(M_{p+1}\)

Pour \(q<p\), on veut tester

\[H_0\; :\; \text{modèle}\; M_{q+1}:\]

\[Y_i=\beta+ \sum_{j=1}^q\alpha_j x_{j,i}+E_i,\quad E_i\; \overset{i.i.d.}{\sim} {\cal N}(0,\sigma^2)\]

\[H_1\; :\; \text{modèle}\; M_{p+1}:\]

\[Y_i= \beta + \sum_{j=1}^p\alpha_j x_{j,i}+E_i,\quad E_i\; \overset{i.i.d.}{\sim} {\cal N}(0,\sigma^2)\]

\[\begin{align}T_n &= \frac{\left(SCR(M_{q+1}) - SCR(M_{p+1})\right)\Bigl/(p-q)}{SCR(M_{p+1})\Bigl/~(n-p-1)}\\& \sim_{H_0}^{}~\mathcal{F}(p-q\,,\,n-p-1)\end{align}\]

5.3 Dans le modèle linéaire simple

Rappel

Dans le modèle \(M_2\) : \(\forall i \in \{1, \dots, n\} \;\;\) \[ Y_i =\beta + \alpha x_i + E_i,\quad E_i \overset{i.i.d.}{\sim}\mathcal{N}(0, \sigma^2)\]

L’estimateur de la variance résiduelle était :

\[\begin{align} S^2 &= \frac{1}{n-2}\sum_{i=1}^n(Y_i-\widehat{Y_i})^2\\ &=\frac{1}{n-2}\sum_{i=1}^n(Y_i-Ax_i-B)^2\\ &= \frac{SCR(M_2)}{n-2}, \end{align}\]

où \(\widehat{Y_i}=Ax_i+B\), la prévision (aléatoire) par le modèle de régression linéaire associée à \(x_i\).

5.4 Le paramètre de variance résiduelle

Objectif

Dans le modèle \(M_{p+1}\):

\(\forall i \in \{1, \dots, n\} \;\;\) \[ Y_i =\beta + \sum_{j=1}^p\alpha_j x_{j,i} + E_i, \quad E_i \overset{i.i.d.}{\sim}\mathcal{N}(0, \sigma^2)\]

il y’a un paramètre de variance \(\sigma ^2\) : la variance residuelle.

Quel est son estimateur ?

5.5 Estimation de la variance residuelle du modèle

Définition : Prévision aléatoire par le modèle \(M_{p+1}\)

\[\begin{align}\widehat{Y_i}&=B+A_1x_{i,1}+\cdots+A_px_{i,p}\\&=B+\sum_{j=1}^p A_jx_{i,j}\end{align}\]

Définition : Résidus aléatoires du modèle \(M_{p+1}\)

\[E_i=Y_i-\widehat{Y_i}\]

Définition : Somme des carrés résiduels

\[\begin{align}SCR(M_{p+1})&=\sum_{i=1}^n E_i^2 \\ &=\sum_{i=1}^n\left(Y_i-(B+\sum_{j=1}^p A_jx_{i,j})\right)^2\end{align}\]

5.6 D’un point de vue matriciel

Définition : Vecteur (aléatoire) des prévisions

\[\widehat{Y}\ =\ \begin{bmatrix} \widehat{Y_1}\\ \widehat{Y_2}\\ \vdots \\ \widehat{Y_n}\end{bmatrix}=\ \begin{bmatrix} B+\sum_{j=1}^p A_jx_{1,j}\\ B+\sum_{j=1}^p A_jx_{2,j}\\ \vdots \\ B+\sum_{j=1}^p A_jx_{n,j}\end{bmatrix}=X\hat{\theta}\]

Définition : Vecteur (aléatoire) des résidus

\[E\ =\ \begin{bmatrix} E_1\\ E_2\\ \vdots \\ E_n\end{bmatrix}=\ \begin{bmatrix} Y_1-\widehat{Y_1}\\ Y_2-\widehat{Y_2}\\ \vdots \\ Y_n-\widehat{Y_n}\end{bmatrix}=Y-X\hat{\theta}\]

Définition : Somme des carrés résiduels

\[\begin{align}SCR(M_{p+1})&=\sum_{i=1}^n E_i^2 \\ &= ||E||^2\\ &=||Y-X\hat{\theta}||^2\end{align}\]

5.7 Estimateur de la variance résiduelle \(\sigma^2\)

Définition : Estimateur de \(\sigma^2\)

\[\begin{align} {S^2}_{(M_{p+1})}&=\frac{SCR(M_{p+1})}{n-(p+1)}\\ &=\frac{1}{n-p-1}\sum_{i=1}^n E_i^2\end{align}\]

\[\begin{align}{S^2}_{(M_{p+1})}&=\frac{\displaystyle\sum_{i=1}^n\left(Y_i-(B+\sum_{j=1}^p A_jx_{i,j})\right)^2}{n-p-1}\\ &=\frac{||Y - X\,\hat{\theta}||^2}{n-p-1} \end{align}\]

Définition : Estimation de \(\sigma^2\)

Réalisation \(s^2\) de \({S^2}_{(M_{p+1})}\) sur les données

\[\begin{align}s^2&=\frac{1}{n-p-1}\sum_{i=1}^n e_i(M_{p+1})^2\\&=\sum_{i=1}^n(y_i-(b+a_1x_{i,1}+\cdots+a_px_{i,p}))^2\end{align}\] où \(e_i(M_{p+1})\) sont les résidus observés dans le modèle \(M_{p+1}\)

6 Propriétés et lois des estimateurs

6.1 Propriétés et lois de \(S^2\)

Théorème

Sous les hypothèses du modèle linéaire gaussien, \({S^2}_{(M_{p+1})}\) est un estimateur sans biais de \(\sigma^2\) et on a \[\begin{align}&\frac{(n-p-1){S^2}_{(M_{p+1})}}{\sigma^2}\\ &=\frac{\sum_{i=1}^n(Y_i-B-\sum_{j=1}^pA_jx_{i,j})^2}{\sigma^2}\\&\sim\chi^2(n-p-1)\end{align}\]

De plus \({S^2}_{(M_{p+1})}\) est indépendant de \(\hat{\theta}\).

6.2 Propriétés et loi de \(\hat{\theta}\)

Théorème

\(\hat{\theta}\) est un vecteur gaussien d’espérance \(\theta\) et de matrice de variance-covariance

\[Var(\hat{\theta})=\sigma^2(X^t\,X)^{-1}\]

estimée par

\[S^2(\hat{\theta})={S^2}_{(M_{p+1})}(X^t\,X)^{-1}\]

On montre que \[\widehat{\theta}\sim{\cal N}(\theta\,,\,\sigma^2(X^t\,X)^{-1})\]

En particulier, pour tout \(k=1,\cdots,p+1\), si \(c_{kk}={(X^t\,X)^{-1}}_{kk}\) désigne le \(k\)-ème élément diagonal de \((X^t\,X)^{-1}\), et \(\widehat{\theta}_k\) désigne le \(k\)-ème élément du vecteur \(\hat{\theta}\), on a \[\widehat{\theta}_k\sim{\cal N}(\theta_k\,,\,\sigma^2c_{kk})\]

En remplaçant \(\sigma^2\) par son estimateur on a, \[\frac{(\widehat{\theta}_k-\theta_k)}{\sqrt{{S^2}_{(M_{p+1})}c_{kk}}}\sim \mathcal{T}{(n-p-1)}\]

6.3 Démonstrations

On rappelle que \[ \hat{\theta}=(X^t\,X)^{-1} X^t Y \]

\[ \begin{align}E(\hat{\theta})&=E((X^t\,X)^{-1} X^t Y)\\ &=(X^t\,X)^{-1} X^t E(Y)\\ &=(X^t\,X)^{-1} X^t X\theta\\ &=\theta \end{align}\]

\[\begin{align} {\rm{Var}}(\hat{\theta})&={\rm{Var}}((X^t\,X)^{-1} X^t Y)\\ &= (X^t\,X)^{-1} X^t {\rm{Var}}(Y) (X^t\,X)^{-1} X^t)^t\\ &= (X^t\,X)^{-1} X^t \sigma^2I_n ((X^t\,X)^{-1} X^t)^t\\ &= \sigma^2(X^t\,X)^{-1} X^t (X^t\,X)^{-1} X^t)^t\\ &= \sigma^2(X^t\,X)^{-1} X^t ((X^t)^t) ((X^t\,X)^{-1})^t \\ &= \sigma^2(X^t\,X)^{-1} X^t X ((X^t\,X)^{t})^{-1} \\ &= \sigma^2 (X^t\,(X^{t})^t)^{-1} \\ &= \sigma^2 (X^t\,X)^{-1} \\ \end{align}\]

6.4 Test et intervalle de confiance

À partir du résultat pour tout \(k=1,\cdots,p+1\)

\[\frac{(\hat{\theta}_k-\theta_k)}{\sqrt{{S^2}_{(M_{p+1})}c_{kk}}}\sim \mathcal{T}{(n-p-1)}\]

On peut construire un intervalle de confiance de chaque paramètre de régression \(\theta_k\), où \(\theta_k\) désigne l’un des paramètres d’espérance (\(\beta\), ou \(\alpha_1\), …, ou \(\alpha_p\)).

On peut construire le test de Student de la nullité de chaque coefficient de régression \(H_0:\theta_k=0\) contre \(H_1:\theta_k\neq 0\) à partir de la statistique de test \[T_n=\frac{\hat{\theta}_k}{\sqrt{{S^2}_{(M_{p+1})}c_{kk}}}\sim_{H_0} \mathcal{T}{(n-p-1)}\]

6.5 Intervalle de confiance de \(\theta_k\)

Théorème : Estimateurs par intervalle de niveau de confiance \(1-\delta\) des \(\theta_k\)

À partir des lois de \(\theta_k\) et \(S^2_{(M_{p+1})}\), on obtient:

\[\begin{align}IC_{1-\delta}(\theta_k) = &\left[\hat{\theta}_k-t_{1-\frac{\delta}{2}} \sqrt{{S^2}_{(M_{p+1})}c_{kk}}\right.;\\ &\left.\;\hat{\theta}_k+t_{1-\frac{\delta}{2}} \sqrt{{S^2}_{(M_{p+1})}c_{kk}}\;\right]\\ \end{align}\]

où \(t_{1-\frac{\delta}{2}}\) est tel que \(\mathbb{P}\left(\mid \mathcal{T}(n-p-1)\mid \leq t_{1-\frac{\delta}{2}}\right)=1-\delta\).

\(t_{1-\frac{\delta}{2}}\) est le quantile d’ordre \(1-\frac{\delta}{2}\) de la loi \(\mathcal{T}(n-p-1)\).

6.6 Intervalle de confiance de \(\theta_k\) avec R :

confint(mod, level = 0.95)

                   2.5 %     97.5 %
(Intercept) -1446.805131 -650.33513
Length1        -9.153972   47.07131
Height        -11.528580  108.52149
Width         -42.941396  176.37258

7 Test dans le modèle de régression linéaire multiple

7.1 Test dans le modèle de régression multiple

Objectifs

Tester si un coefficient \(\theta_k\) est significativement différent de \(0\).

Tester si la contribution globale des variables explicatives est significative pour expliquer \(Y\).

Déterminer la contribution effective des variables explicatives dans la modélisation de \(Y\).

Tester si un ensemble de \(q\) variables explicatives (\(q<p\)) ne suffit pas à modéliser correctement \(Y\).

7.2 Test de Student pour \(\theta_k\)

Statistique de test et loi sous \(H_0\)

\[{T_n= \frac{\hat{\theta}_k}{\sqrt{{S^2}_{(M_{p+1})}c_{kk}}}\sim_{H_0}\mathcal{T}(n- p -1)}\]

Zone de rejet

\[R_\delta = \{|T|> t_{1-\frac{\delta}{2}} \} \]

On rejette \(H_0\) si \(t_{obs}> t_{1-\frac{\delta}{2}}\)

Application numérique

On calcule \(t_{obs} = \frac{a}{\sqrt{{s^2}_{(M_{p+1})}c_{kk}}}\) la réalisation de \(T_n\).

On compare avec \(t_{1-\frac{\delta}{2}}\) et on conclut.

7.3 \(p\)-valeur du test de Student.

Calcul de la \(p\)-valeur

\[\begin{align}p_c&=\mathbb{P}_{H_0}(\mid T_n\mid >\mid t_n\mid)\\ &=2(1-\mathbb{P}(T_n\leq |t_n|))\end{align}\] où \(T_n\sim \mathcal{T}(n-p -1)\)

Inteprétation

Pour un risque de 1ère espèce \(\delta\) fixé acceptable (par ex \(\delta=5\%\))

si \(p_c <\delta\), le test de niveau \(\delta\) est significatif (liaison significative), on rejette \(H_0\).
si \(p_c >\delta\), le test de niveau \(\delta\) n’est pas significatif (liaison non significative), on ne rejette pas \(H_0\).

7.4 Test de Student avec R

mod <- lm(Weight ~Length1 + Height + Width , data =  fish)
summary(mod)


Call:
lm(formula = Weight ~ Length1 + Height + Width, data = fish)

Residuals:
     Min       1Q   Median       3Q      Max 
-183.872  -17.044   -0.495   24.591  101.032 

Coefficients:
            Estimate Std. Error t value Pr(>|t|)    
(Intercept) -1048.57     187.86  -5.582 4.13e-05 ***
Length1        18.96      13.26   1.430    0.172    
Height         48.50      28.31   1.713    0.106    
Width          66.72      51.73   1.290    0.215    
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Residual standard error: 63.03 on 16 degrees of freedom
Multiple R-squared:  0.808, Adjusted R-squared:  0.772 
F-statistic: 22.44 on 3 and 16 DF,  p-value: 5.624e-06

7.5 Test de la contribution globale des variables explicatives

Généralisation du test de Fisher à \(p\) variables explicatives

\(H_0\) : aucune des \(p\) variables explicatives n’a d’influence sur \(Y\).

\(H_1\) : au moins une des variables explicatives contribue à expliquer \(Y\).

Autre formulation :

Test du modèle constant \[H_0\ :\ \text{modèle}\ M_1:\]

\[Y_i=\beta+E_i,\quad E_i\ \overset{i.i.d.}{\sim}\ {\cal N}(0,\sigma^2)\] contre le modèle complet

\[H_1\ :\ \text{modèle}\ M_{p+1}:\]

\[Y_i=\beta+\sum_{j=1}^p\alpha_j x_{i,j}+E_i,\quad E_i\ \overset{i.i.d.}{\sim} {\cal N}(0,\sigma^2)\]

Autre formulation :

Test de \(H_0:\ \forall j=1,\cdots,p,\ \alpha_j=0\quad\) contre \(\quad H_1:\exists j\ \alpha_j\neq 0\)

7.6 Test de Fisher de la contribution globale des variables explicatives

Définition : \(SCT\)

La variabilité de \(Y\) sans tenir compte du modèle.

\[\color{purple}{SCT =\displaystyle\sum_{i = 1}^n( Y_i - \bar{Y})^2}\]

Définition : \(SCM(M_{p+1})\)

Partie de la variabilité de \(Y\) expliquée par le modèle \(M_{p+1}\).

\[\color{blue}{SCM(M_{p+1}) = \displaystyle\sum_{i=1}^n(\widehat{Y_i}(M_{p+1})-\bar{Y})^2}\]

Définition : \(SCR(M_{p+1})\)

Partie de la variabilité de \(Y\) qui n’est pas expliquée par le modèle \(M_{p+1}\).

\[\color{red}{\begin{align}SCR(M_{p+1}) &= \displaystyle\sum_{i=1}^n(Y_i-\widehat{Y_i}(M_{p+1}))^2\\&=\displaystyle\sum_{i=1}^n E_i(M_{p+1}) ^2\end{align}}\]

Théorème de décomposition de la variance du modèle \(M_{p+1}\)

\[\color{purple}{SCT}=\color{blue}{SCM(M_{p+1})}+ \color{red}{SCR(M_{p+1})}\]

7.7 Test de Fisher de la contribution globale des variables explicatives (suite)

Notations

Pour mieux distinguer \(M_1\) et \(M_{p+1}\), nous allons les noter \(\textcolor{green}{M_1}\) et \(\textcolor{magenta}{M_{p+1}}\).

Remarques

Dans le modèle \(\textcolor{green}{M_1}\), \(\beta\) est estimée par \(\textcolor{purple}{\bar{Y}}\) et \(\widehat{Y_i}(\textcolor{green}{M_{1}})=\textcolor{purple}{\bar{Y}}\) donc

\[\begin{align}SCR(\textcolor{green}{M_1})&=\sum_{i=1}^n(Y_i-\widehat{Y_i}(\textcolor{green}{M_{1}}))^2\\ &=\textcolor{purple}{\sum_{i=1}^n(Y_i-\bar{Y})^2}\\ &=\textcolor{purple}{SCT}\end{align}\]

On a donc \(SCM(\textcolor{magenta}{M_{p+1}})=SCR(\textcolor{green}{M_1})-SCR(\textcolor{magenta}{M_{p+1}})\) : représente la réduction d’erreur quand on passe du modèle \(\textcolor{green}{M_1}\) au modèle \(\textcolor{magenta}{M_{p+1}}\).

7.8 Test de \(H_0\) : modèle \(M_1\) contre \(H_1\) : modèle \(M_{p+1}\)

Statistique de test

\[\begin{align}T_n&=\frac{SCM(\textcolor{magenta}{M_{p+1}})/p}{SCR(\textcolor{magenta}{M_{p+1}})/(n-p-1)}\\&\sim_{H_0} \mathcal{F}(p,n-p-1)\end{align}\]

qui s’écrit donc aussi

\[\begin{align}T_n&=\frac{(SCR(\textcolor{green}{M_1})-SCR(\textcolor{magenta}{M_{p+1}}))/p}{SCR(\textcolor{magenta}{M_{p+1}})/(n-p-1)}\\&\sim_{H_0} \mathcal{F}(p,n-p-1)\end{align}\]

Zone de rejet

\[R_\delta = \{T_n > f_{1-\delta} \}\]

\(f_{1-\delta}\) est le quantile \(1 - \delta\) de la loi de Fisher \(\mathcal{F}(p,n-p-1)\).

Si \(t_{obs} > f_{1-\delta}\), on rejette \(H_0\) au risque \(\delta\).

7.9 Test de \(H_0\) : modèle \(M_1\) contre \(H_1\) : modèle \(M_{p+1}\)

Mise en oeuvre et conclusion

Si \(t_{obs} > f_{1-\delta}\), on conserve le modèle complet \(\textcolor{magenta}{M_{p+1}}\) et la contribution globale des \(p\) variables explicatives est significative : au moins une des variables explicatives contribue à expliquer la variabilité de \(Y\).
Si \(t_{obs} < f_{1-\delta}\), on ne rejette pas \(H_0\) : le modèle \(\textcolor{green}{M_1}\) suffit et les variables explicatives ne contribuent pas à expliquer de manière significative \(Y\).

\(p\)-valeur

\[p_c=\mathbb{P}_{H_0}(T_n > t_{obs})\] que l’on compare à un risque, par exemple \(\delta=5\%\)

On rejette \(H_0\) si \(p_c<\delta\).

7.10 Test avec R

mod <- lm(Weight ~Length1 + Height + Width , data =  fish)
summary(mod)


Call:
lm(formula = Weight ~ Length1 + Height + Width, data = fish)

Residuals:
     Min       1Q   Median       3Q      Max 
-183.872  -17.044   -0.495   24.591  101.032 

Coefficients:
            Estimate Std. Error t value Pr(>|t|)    
(Intercept) -1048.57     187.86  -5.582 4.13e-05 ***
Length1        18.96      13.26   1.430    0.172    
Height         48.50      28.31   1.713    0.106    
Width          66.72      51.73   1.290    0.215    
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Residual standard error: 63.03 on 16 degrees of freedom
Multiple R-squared:  0.808, Adjusted R-squared:  0.772 
F-statistic: 22.44 on 3 and 16 DF,  p-value: 5.624e-06

7.11 Table de l’analyse de la variance de la régression

\[M_1:\quad Y_i=\beta+ E_i, \quad E_i \overset{i.i.d.}{\sim} {\cal N}(0\,,\,\sigma^2)\]

m1 <- lm(Weight ~1, data = fish)
anova(m1, mod)

Analysis of Variance Table

Model 1: Weight ~ 1
Model 2: Weight ~ Length1 + Height + Width
  Res.Df    RSS Df Sum of Sq      F    Pr(>F)    
1     19 331013                                  
2     16  63568  3    267445 22.439 5.624e-06 ***
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Comme le test global \(H_0\) : \(M_1\) contre \(H_1\) : \(M_{p+1}\) est significatif, au moins une des variables explicatives contribue à expliquer \(Y\).

7.12 Test du modèle \(M_{q+1}\) contre le modèle \(M_{p+1}\)

Objectif

Tester si un ensemble de \(\textcolor{green}{q}\) variables explicatives ne suffit pas à expliquer \(Y\), avec \(\textcolor{green}{q}<\textcolor{magenta}{p}\).

Autre formulation :

\[H_0\ :\ \text{modèle}\ \textcolor{green}{M_{q+1}}\ :\]

\[Y_i=\beta+\sum_{j=1}^{\textcolor{green}{q}}\alpha_j x_{i,j}+E_i,\quad E_i\overset{i.i.d.}{\sim}{\cal N}(0,\sigma^2)\] contre l’alternative \[H_1\ :\ \text{modèle}\ \textcolor{magenta}{M_{p+1}}:\]

\[Y_i=\beta+\sum_{j=1}^{\textcolor{magenta}{p}}\alpha_j x_{i,j}+E_i,\quad E_i\ \overset{i.i.d.}{\sim} {\cal N}(0,\sigma^2)\]

Autre formulation :

Test de \(H_0:\ \forall j=\textcolor{green}{q}+1,\cdots,\textcolor{magenta}{p},\ \alpha_j=0 \quad\) contre \(\quad H_1:\exists j=\textcolor{green}{q}+1,\cdots,\textcolor{magenta}{p}\ \alpha_j\neq 0\)

7.13 Description des deux modèles \(M_{q+1}\) et \(M_{p+1}\)

Modèle à \(\textcolor{magenta}{p}\) variables explicatives \(\textcolor{magenta}{M_{p+1}}\)

\[ Y_i = \beta\ +\ \sum_{j=1}^{\textcolor{magenta}{p}}\alpha_j\,x_{i,j}+ E_i,\] \[E_i\ \overset{i.i.d.}{\sim}\ {\cal N}(0,\sigma^2).\]

Sous la forme matricielle :

\[(\textcolor{magenta}{M_{p+1}})\,: Y = X^{(\textcolor{magenta}{p})}\theta^{(\textcolor{magenta}{p})}\ +\ E,\]

\[E \sim {\cal N}(0, \sigma^2 I_n) \]

avec

\[ X^{(\textcolor{magenta}{p})} = \begin{bmatrix} 1 & x_{1,1} & \ldots & x_{1,\textcolor{magenta}{p}}\\ 1 & x_{2,1} & \ldots & x_{2,\textcolor{magenta}{p}}\\ \vdots & \vdots & & \vdots \\ 1 & x_{n,1} & \ldots & x_{n,\textcolor{magenta}{p}} \end{bmatrix} \]

\[\theta^{(\textcolor{magenta}{p})}=[\beta,\alpha_1,\cdots,\alpha_\textcolor{magenta}{p}]^t\]

7.14 Description des deux modèles \(M_{q+1}\) et \(M_{p+1}\)

Modèle à \(\textcolor{green}{q}\) variables explicatives \(\textcolor{green}{M_{q+1}}\)

\[ Y_i = \beta\ +\ \sum_{j=1}^{\textcolor{green}{q}}\alpha_j\,x_{i,j}+ E_i\]

\[E_i\ \overset{i.i.d.}{\sim}\ {\cal N}(0,\sigma^2).\]

Sous la forme matricielle :

\[(\textcolor{green}{M_{q+1}})\,: Y = X^{(\textcolor{green}{q})}\theta^{(\textcolor{green}{q})}\ +\ E\]

\[E \sim {\cal N}(0, \sigma^2 I_n) \]

avec

\[ X^{(\textcolor{green}{q})} = \begin{bmatrix} 1 & x_{1,1} & \ldots & x_{1,\textcolor{green}{q}}\\ 1 & x_{2,1} & \ldots & x_{2,\textcolor{green}{q}}\\ \vdots & \vdots & & \vdots \\ 1 & x_{n,1} & \ldots & x_{n,\textcolor{green}{q}} \end{bmatrix} \] et \[\theta^{(\textcolor{green}{q})}=[\beta,\alpha_1,\cdots,\alpha_\textcolor{green}{q}]^t\]

Le modèle à \(\textcolor{green}{q}\) variables explicatives \(\textcolor{green}{M_{q+1}}\) est emboîté dans le modèle à \(\textcolor{magenta}{p}\) variables explicatives \(\textcolor{magenta}{M_{p+1}}\) : les \(\textcolor{green}{q}\) variables explicatives forment un sous-ensemble des \(\textcolor{magenta}{p}\) variables explicatives.

7.15 Estimateur des paramètres dans \(M_{p+1}\) et dans \(M_{q+1}\)

Paramètres d’espérances

\[\begin{align} \hat{\theta}^{(\textcolor{magenta}{p})} &= (B^{(\textcolor{magenta}{p})},A_1^{(\textcolor{magenta}{p})},\dots,A_{\textcolor{magenta}{p}}^{(\textcolor{magenta}{p})})^t\\ &= \left((X^{(\textcolor{magenta}{p})})^t\,X^{(\textcolor{magenta}{p})}\right)^{-1}\,\!(X^{(\textcolor{magenta}{p})})^t\,Y\\ \hat{\theta}^{(\textcolor{green}{q})} & = (B^{(\textcolor{green}{q})},A_1^{(\textcolor{green}{q})},\dots,A_{\textcolor{green}{q}}^{(\textcolor{green}{q})})^t\\ &=\left((X^{(\textcolor{green}{q})})^t\,X^{(\textcolor{green}{q})}\right)^{-1}\,\!(X^{(\textcolor{green}{q})})^t\,Y\ \end{align}\]

Remarques

\(B^{(\textcolor{magenta}{p})}\neq B^{(\textcolor{green}{q})}, A_1^{(\textcolor{magenta}{p})}\neq A_1^{(\textcolor{green}{q})},\cdots\)

7.16 Prévision et résidu dans \(M_{p+1}\) et dans \(M_{q+1}\)

Prévision dans \(M_{p+1}\) et dans \(M_{q+1}\)

\[\begin{align}\widehat{Y_i}(\textcolor{magenta}{M_{p+1}})&=B^{(\textcolor{magenta}{p})}+\sum_{j=1}^{\textcolor{magenta}{p}} A^{(\textcolor{magenta}{p})}_jx_{i,j}\\ \widehat{Y_i}(\textcolor{green}{M_{q+1}})&=B^{(\textcolor{green}{q})}+\sum_{j=1}^{\textcolor{green}{q}} A^{(\textcolor{green}{q})}_jx_{i,j} \end{align}\]

Résidu dans \(M_{p+1}\) et dans \(M_{q+1}\)

\[\begin{align}E_i(\textcolor{magenta}{M_{p+1}})&=Y_i-\widehat{Y_i}(\textcolor{magenta}{M_{p+1}})\\ &=Y_i-(B^{(\textcolor{magenta}{p})}+\sum_{j=1}^{\textcolor{magenta}{p}} A^{(\textcolor{magenta}{p})}_jx_{i,j})\\ E_i(\textcolor{green}{M_{q+1}})&=Y_i-\widehat{Y_i}(\textcolor{green}{M_{q+1}})\\ &=Y_i-(B^{(\textcolor{green}{q})}+\sum_{j=1}^{\textcolor{green}{q}} A^{(\textcolor{green}{q})}_jx_{i,j}) \end{align}\]

7.17 SCR dans \(M_{p+1}\) et dans \(M_{q+1}\)

SCR dans \(M_{p+1}\) et dans \(M_{q+1}\)

\[\begin{align}SCR(\textcolor{magenta}{M_{p+1}}) &= \sum_{i=1}^{n}\left(Y_i-B^{(\textcolor{magenta}{p})}-\sum_{j=1}^{p}A^{(\textcolor{magenta}{p})}_j\,x_{i,j}\right)^2\\ &=\ ||Y - X^{(\textcolor{magenta}{p})}\hat{\theta}^{(\textcolor{magenta}{p})}||^2\end{align}\]

\[\begin{align}SCR(\textcolor{green}{M_{q+1}})& = \sum_{i=1}^{n}\left(Y_i-B^{(\textcolor{green}{q})}-\sum_{j=1}^{\textcolor{green}{q}}A^{(\textcolor{green}{q})}_j\,x_{i,j}\right)^2\\& =\ ||Y - X^{(\textcolor{green}{q})}\hat{\theta}^{(\textcolor{green}{q})}||^2\end{align}\]

Estimation de \(\sigma^2\) dans \(M_{p+1}\) et dans \(M_{q+1}\)

\[S^2_{(\textcolor{magenta}{M_{p+1}})}=\frac{SCR(\textcolor{magenta}{M_{p+1}})}{n\textcolor{magenta}{-p-1}}\]

\[S^2_{(\textcolor{green}{M_{q+1}})}=\frac{SCR(\textcolor{green}{M_{q+1}})}{n\textcolor{green}{-q-1}}\]

7.18 Résultat important

Théorème important du cours

Sous les hypothèses du modèle linéaire Gaussien multiple, et si on définit \(SCE=SCR(\textcolor{green}{M_{q+1}})-SCR(\textcolor{magenta}{M_{p+1}})\) alors

Sous le modèle \(\textcolor{green}{M_{q+1}}\), c’est à dire sous \(H_0\) :

\[\frac{SCR(\textcolor{green}{M_{q+1}})}{\sigma^2} \underset{H_0}{\sim}\,\chi^2_{n\textcolor{green}{-q-1}}\]

\[\frac{SCE}{\sigma^2} \underset{H_0}{\sim}\,\chi^2_{\textcolor{magenta}{p}-\textcolor{green}{q}}\\\]

et \(SCE\) et \(SCR(\textcolor{magenta}{M_{p+1}})\) sont indépendants ce qui permet de déduire que

\[\begin{align}&\frac{\left(SCR(\textcolor{green}{M_{q+1}}) - SCR(\textcolor{magenta}{M_{p+1}})\right)\Bigl/(\textcolor{magenta}{p}-\textcolor{green}{q})}{SCR(\textcolor{magenta}{M_{p+1}})\Bigl/~(n\textcolor{magenta}{-p-1})}\\ &\sim_{H_0}^{}~\mathcal{F}(\textcolor{magenta}{p}-\textcolor{green}{q}\,,\,n\textcolor{magenta}{-p-1})\end{align}\]

7.19 Test de \(H_0\) : modèle \(M_{q+1}\) contre \(H_1\) : modèle \(M_{p+1}\)

Statistique de test

\[\begin{align}T_n&=\frac{(SCR(\textcolor{green}{M_{q+1}})-SCR(\textcolor{magenta}{M_{p+1}}))/(\textcolor{green}{p}-\textcolor{magenta}{q})}{SCR(\textcolor{magenta}{M_{p+1}})/(n\textcolor{magenta}{-p-1})}\\ &\sim_{H_0} \mathcal{F}(\textcolor{magenta}{p}-\textcolor{green}{q},n\textcolor{magenta}{-p-1})\end{align}\]

\(SCR(\textcolor{green}{M_{q+1}})-SCR(\textcolor{magenta}{M_{p+1}})\) représente la réduction d’erreur quand on passe de \(\textcolor{green}{M_{q+1}}\) à \(\textcolor{magenta}{M_{p+1}}\).

Zone de rejet

\[R_\delta = \{T_n > f_{1-\delta} \} \] où \(f_\delta\) est le quantile d’ordre \(1-\delta\) de la \(\mathcal{F}(p-q,n-p-1)\).

On rejette \(H_0\) si \(t_{n}> f_{1-\delta}\)

7.20 Test de \(H_0 : \mbox{ modèle } M_{q+1}\) contre \(H_1\) : modèle \(M_{p+1}\)

Interprétation

Si \(t_n > f_{1-\delta}\), on conserve le modèle complet \(\textcolor{magenta}{M_{p+1}}\), on considère que le passage de \(\textcolor{green}{M_{q+1}}\) à \(\textcolor{magenta}{M_{p+1}}\) est significatif : au moins une variable explicative parmi \(x_{\textcolor{green}{q}+1},\cdots,x_{\textcolor{magenta}{p}}\) a une influence significative sur \(Y\) (en plus de \(x_1,\cdots,x_q\)).

Si \(t_n \leq f_{1-\delta}\), on ne rejette pas \(H_0\), on considère que le passage de \(\textcolor{green}{M_{q+1}}\) à \(\textcolor{magenta}{M_{p+1}}\) n’est pas significatif : l’influence des variables explicatives \(x_{\textcolor{green}{q}+1},\cdots,x_{\textcolor{magenta}{p}}\) n’est pas significative pour expliquer \(Y\).

Calcul de la \(p\)-valeur

\[\begin{align}p_c&=\mathbb{P}_{H_0}(T_n > t_n)\\&=\mathbb{P}(\mathcal{F}(p-q,n-p-1)>t_n)\end{align}\]

7.21 Table de l’analyse de la variance de la régression

Avec R, commande anova

\[\begin{array}{|c|c|c|c|c|c|c|} \hline Source & ddl & SCR& ddl & SCE\quad expliquée& Statistique & p-value\\ & résiduelle& && par\quad \textcolor{magenta}{M_{p+1}} & de\quad test & \\ \hline \text{Modèle} & n\textcolor{green}{-q-1} & SCR(\textcolor{green}{M_{q+1}}) & &&&\\ \textcolor{green}{M_{q+1}} &&&&&&\\ \hline \text{Modèle} & n\textcolor{magenta}{-p-1} & SCR(\textcolor{magenta}{M_{p+1}}) & \textcolor{green}{p}-\textcolor{magenta}{q} & SCE=SCR(\textcolor{green}{M_{q+1}}) & t_n=&\mathbb{P}(\mathcal{F}(\textcolor{magenta}{p}-\textcolor{green}{q},n\textcolor{magenta}{-p-1})>t_n)\\ \textcolor{magenta}{M_{p+1}} &&&& - SCR(\textcolor{magenta}{M_{p+1}}) &\frac{SCE/(\textcolor{magenta}{p}-\textcolor{green}{q})}{SCR(\textcolor{magenta}{M_{p+1}})/(n\textcolor{magenta}{-p-1})} &\\ \hline \end{array}\]

7.22 Retour à l’exemple

On va tester

\[\begin{align} &\textcolor{magenta}{M_4} : Y_i=\beta+ \alpha_1x_{i,1}+\alpha_2x_{i,2}+\alpha_3x_{i,3}+ E_i \\ &\textcolor{green}{M_3} : Y_i=\beta+ \alpha_1x_{i,1}+\alpha_2x_{i,2} +E_i \end{align}\]

On va définir le modèle \(\textcolor{green}{M_3}\)

m3 <- lm(Weight ~ Length1 + Width, data = fish)
summary(m3)


Call:
lm(formula = Weight ~ Length1 + Width, data = fish)

Residuals:
     Min       1Q   Median       3Q      Max 
-211.090  -14.742    2.901   25.347  104.568 

Coefficients:
            Estimate Std. Error t value Pr(>|t|)    
(Intercept) -998.742    195.859  -5.099 8.91e-05 ***
Length1       35.713      9.449   3.779   0.0015 ** 
Width         97.835     51.111   1.914   0.0726 .  
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Residual standard error: 66.52 on 17 degrees of freedom
Multiple R-squared:  0.7728,    Adjusted R-squared:  0.746 
F-statistic:  28.9 on 2 and 17 DF,  p-value: 3.391e-06

7.23 Sorties R

anova(m3,  mod)

Analysis of Variance Table

Model 1: Weight ~ Length1 + Width
Model 2: Weight ~ Length1 + Height + Width
  Res.Df   RSS Df Sum of Sq      F Pr(>F)
1     17 75223                           
2     16 63568  1     11655 2.9335 0.1061

anova(m1,  m3)

Analysis of Variance Table

Model 1: Weight ~ 1
Model 2: Weight ~ Length1 + Width
  Res.Df    RSS Df Sum of Sq      F    Pr(>F)    
1     19 331013                                  
2     17  75223  2    255790 28.904 3.391e-06 ***
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

7.24 Remarques sur les modèles emboités

Remarques

On peut toujours comparer des modèles emboités l’un dans l’autre.
On peut donc comparer tous les modèles au modèle complet (structure la plus riche).
On peut de même comparer tous les modèles au modèle constant (structure la moins riche).
On ne peut pas comparer des modèles non emboîtés à l’aide de ce test.

7.25 Ce qu’il faut retenir de ce cours 1/2

Estimateur de \(\sigma^2\) et loi

\[{S^2}_{(M_{p+1})}=\frac{SCR(M_{p+1})}{n-(p+1)}=\frac{||Y - X\,\widehat{\theta}||^2}{n-p-1}\]

\[\frac{(n-p-1){S^2}_{(M_{p+1})}}{\sigma^2}\sim\chi^2(n-p-1)\]

Loi de \(\widehat{\theta}\) et \(\widehat{\theta}_k\)

\[\widehat{\theta} \sim \mathcal{N}\left(\theta,\sigma^2(X^tX)^{-1}\right)\] \[\widehat{\theta}_k \sim \mathcal{N}\left(\theta_k,\sigma^2c_{kk}\right)\text{ avec } c_{kk} = (X^tX)^{-1}_{k,k}\]

\[\frac{\widehat{\theta}_k- \theta_k}{\sqrt{S^2_{(M_{p+1})}c_{kk}}}\sim\mathcal{T}(n-p-1) \]

Intervalle de confiance de \(\theta_k\)

7.26 Ce qu’il faut retenir de ce cours 2/2

Test de Student de la nullité de la pente

\[H_0:\alpha_k=0\]

\[H_1:\alpha_k\neq 0\]

Test de Fisher de comparaison de \(M_1\) contre \(M_{p+1}\)

\[H_0\; :\; \text{modèle}\; M_1:\] \[Y_i=\beta+E_i,\quad E_i\; \overset{i.i.d.}{\sim} {\cal N}(0,\sigma^2)\]

\[H_1\; :\; \text{modèle}\; M_{p+1}:\] \[Y_i= \beta + \sum_{j=1}^p\alpha_j x_{j,i}+E_i,\; E_i\; \overset{i.i.d.}{\sim} {\cal N}(0,\sigma^2)\]

Test de Fisher de comparaison de \(M_{q+1}\) contre \(M_{p+1}\)

Pour \(q<p\), on veut tester

\[H_0\; :\; \text{modèle}\; M_{q+1}:\]

\[Y_i=\beta+ \sum_{j=1}^q\alpha_j x_{j,i}+E_i,\quad E_i\; \overset{i.i.d.}{\sim} {\cal N}(0,\sigma^2)\]

\[H_1\; :\; \text{modèle}\; M_{p+1}:\]

\[Y_i= \beta + \sum_{j=1}^p\alpha_j x_{j,i}+E_i,\quad E_i\; \overset{i.i.d.}{\sim} {\cal N}(0,\sigma^2)\]

\[\begin{align}T_n &= \frac{\left(SCR(M_{q+1}) - SCR(M_{p+1})\right)\Bigl/(p-q)}{SCR(M_{p+1})\Bigl/~(n-p-1)}\\& \sim_{H_0}^{}~\mathcal{F}(p-q\,,\,n-p-1)\end{align}\]